# from google.colab import drive
# drive.mount('/content/drive')
!pip install -e .
Obtaining file:///Users/M286333/Documents/_projects/hifuku
Preparing metadata (setup.py) ... done
Requirement already satisfied: opencv-python in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from hifuku==0.0.1) (4.6.0)
Requirement already satisfied: albumentations in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from hifuku==0.0.1) (1.3.0)
Requirement already satisfied: pytorch-lightning in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from hifuku==0.0.1) (2.0.0)
Requirement already satisfied: segmentation-models-pytorch in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from hifuku==0.0.1) (0.3.1)
Requirement already satisfied: scikit-image in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from hifuku==0.0.1) (0.19.3)
Requirement already satisfied: opencv-python-headless>=4.1.1 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from albumentations->hifuku==0.0.1) (4.7.0.72)
Requirement already satisfied: qudida>=0.0.4 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from albumentations->hifuku==0.0.1) (0.0.4)
Requirement already satisfied: scipy in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from albumentations->hifuku==0.0.1) (1.9.3)
Requirement already satisfied: numpy>=1.11.1 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from albumentations->hifuku==0.0.1) (1.23.5)
Requirement already satisfied: PyYAML in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from albumentations->hifuku==0.0.1) (6.0)
Requirement already satisfied: imageio>=2.4.1 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from scikit-image->hifuku==0.0.1) (2.23.0)
Requirement already satisfied: PyWavelets>=1.1.1 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from scikit-image->hifuku==0.0.1) (1.4.1)
Requirement already satisfied: tifffile>=2019.7.26 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from scikit-image->hifuku==0.0.1) (2022.10.10)
Requirement already satisfied: pillow!=7.1.0,!=7.1.1,!=8.3.0,>=6.1.0 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from scikit-image->hifuku==0.0.1) (9.4.0)
Requirement already satisfied: networkx>=2.2 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from scikit-image->hifuku==0.0.1) (2.8.8)
Requirement already satisfied: packaging>=20.0 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from scikit-image->hifuku==0.0.1) (22.0)
Requirement already satisfied: fsspec[http]>2021.06.0 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from pytorch-lightning->hifuku==0.0.1) (2022.11.0)
Requirement already satisfied: torchmetrics>=0.7.0 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from pytorch-lightning->hifuku==0.0.1) (0.11.0)
Requirement already satisfied: torch>=1.11.0 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from pytorch-lightning->hifuku==0.0.1) (1.13.0)
Requirement already satisfied: lightning-utilities>=0.7.0 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from pytorch-lightning->hifuku==0.0.1) (0.8.0)
Requirement already satisfied: typing-extensions>=4.0.0 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from pytorch-lightning->hifuku==0.0.1) (4.4.0)
Requirement already satisfied: tqdm>=4.57.0 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from pytorch-lightning->hifuku==0.0.1) (4.64.1)
Requirement already satisfied: timm==0.4.12 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from segmentation-models-pytorch->hifuku==0.0.1) (0.4.12)
Requirement already satisfied: torchvision>=0.5.0 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from segmentation-models-pytorch->hifuku==0.0.1) (0.14.0a0)
Requirement already satisfied: efficientnet-pytorch==0.7.1 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from segmentation-models-pytorch->hifuku==0.0.1) (0.7.1)
Requirement already satisfied: pretrainedmodels==0.7.4 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from segmentation-models-pytorch->hifuku==0.0.1) (0.7.4)
Requirement already satisfied: munch in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from pretrainedmodels==0.7.4->segmentation-models-pytorch->hifuku==0.0.1) (2.5.0)
Requirement already satisfied: requests in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from fsspec[http]>2021.06.0->pytorch-lightning->hifuku==0.0.1) (2.28.2)
Requirement already satisfied: aiohttp!=4.0.0a0,!=4.0.0a1 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from fsspec[http]>2021.06.0->pytorch-lightning->hifuku==0.0.1) (3.8.3)
Requirement already satisfied: scikit-learn>=0.19.1 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from qudida>=0.0.4->albumentations->hifuku==0.0.1) (1.2.0)
Requirement already satisfied: yarl<2.0,>=1.0 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from aiohttp!=4.0.0a0,!=4.0.0a1->fsspec[http]>2021.06.0->pytorch-lightning->hifuku==0.0.1) (1.8.2)
Requirement already satisfied: multidict<7.0,>=4.5 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from aiohttp!=4.0.0a0,!=4.0.0a1->fsspec[http]>2021.06.0->pytorch-lightning->hifuku==0.0.1) (6.0.4)
Requirement already satisfied: async-timeout<5.0,>=4.0.0a3 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from aiohttp!=4.0.0a0,!=4.0.0a1->fsspec[http]>2021.06.0->pytorch-lightning->hifuku==0.0.1) (4.0.2)
Requirement already satisfied: charset-normalizer<3.0,>=2.0 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from aiohttp!=4.0.0a0,!=4.0.0a1->fsspec[http]>2021.06.0->pytorch-lightning->hifuku==0.0.1) (2.1.1)
Requirement already satisfied: attrs>=17.3.0 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from aiohttp!=4.0.0a0,!=4.0.0a1->fsspec[http]>2021.06.0->pytorch-lightning->hifuku==0.0.1) (22.2.0)
Requirement already satisfied: frozenlist>=1.1.1 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from aiohttp!=4.0.0a0,!=4.0.0a1->fsspec[http]>2021.06.0->pytorch-lightning->hifuku==0.0.1) (1.3.3)
Requirement already satisfied: aiosignal>=1.1.2 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from aiohttp!=4.0.0a0,!=4.0.0a1->fsspec[http]>2021.06.0->pytorch-lightning->hifuku==0.0.1) (1.3.1)
Requirement already satisfied: threadpoolctl>=2.0.0 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from scikit-learn>=0.19.1->qudida>=0.0.4->albumentations->hifuku==0.0.1) (3.1.0)
Requirement already satisfied: joblib>=1.1.1 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from scikit-learn>=0.19.1->qudida>=0.0.4->albumentations->hifuku==0.0.1) (1.2.0)
Requirement already satisfied: six in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from munch->pretrainedmodels==0.7.4->segmentation-models-pytorch->hifuku==0.0.1) (1.16.0)
Requirement already satisfied: urllib3<1.27,>=1.21.1 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from requests->fsspec[http]>2021.06.0->pytorch-lightning->hifuku==0.0.1) (1.26.13)
Requirement already satisfied: idna<4,>=2.5 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from requests->fsspec[http]>2021.06.0->pytorch-lightning->hifuku==0.0.1) (3.4)
Requirement already satisfied: certifi>=2017.4.17 in /Users/M286333/miniforge3/envs/py39/lib/python3.9/site-packages (from requests->fsspec[http]>2021.06.0->pytorch-lightning->hifuku==0.0.1) (2023.5.7)
Installing collected packages: hifuku
Attempting uninstall: hifuku
Found existing installation: hifuku 0.0.1
Uninstalling hifuku-0.0.1:
Successfully uninstalled hifuku-0.0.1
Running setup.py develop for hifuku
Successfully installed hifuku-0.0.1
# %cd /content/drive/MyDrive/hifuku
# !pip install -r requirements.txt
import numpy as np
import pandas as pd
import seaborn as sns
import matplotlib.pyplot as plt
import os
from glob import glob
import csv
import os.path
from pathlib import Path
import random
import math
import warnings
warnings.simplefilter('ignore')
import cv2
import PIL
PIL.Image.MAX_IMAGE_PIXELS = 933120000
from PIL import Image
import joblib
from sklearn.manifold import TSNE
from sklearn.linear_model import LinearRegression
from sklearn.preprocessing import PolynomialFeatures
from sklearn.preprocessing import StandardScaler
from sklearn.mixture import GaussianMixture
from sklearn.cluster import KMeans
from scipy import stats
import ast
from scipy.spatial import distance
# from scipy.spatial import Voronoi, voronoi_plot_2d
from scipy.stats import entropy
from sklearn.neighbors import KernelDensity
from shapely.geometry import Polygon
from tqdm import tqdm
import hifuku
palette = sns.color_palette('viridis', n_colors=5)
sns.set(context='talk', style='ticks', palette=palette, font_scale=.7,rc={
'axes.spines.right': False,
'axes.spines.top': False,
'figure.figsize': [4.0, 4.0],
'legend.fancybox': False,
'legend.edgecolor': 'black',
'legend.fancybox': False,
'legend.frameon': True,
'legend.handletextpad': 0.5,
'legend.loc': 'upper right',
'legend.loc': 'best',
'legend.title_fontsize': None,
})
sns.palplot(palette)
colors = [(int(color[0]*255), int(color[1]*255), int(color[2]*255)) for color in palette]
# root = '/content/drive/MyDrive/hifuku'
root = '/Users/M286333/Documents/_projects/hifuku'
# a sample WSI image
is_wsi = True
path = f'{root}/data/sample_wsi_scale_0273.jpg'
scale = 0.2730 # µm / pixel
# if you want to analyze a non-wsi sample image, activate the following 3 lines
# is_wsi = False
# path = f'{root}/data/sample_scale_02325.jpg'
# scale = 0.2325 # µm / pixel
image_id = os.path.splitext(os.path.splitext(os.path.basename(path))[0])[0]
save_dir= f'{root}/results/{image_id}'
image_id
'sample_wsi_scale_0273'
hifuku.main(root, path, is_wsi=is_wsi, scale=scale)
Global seed set to 0
<Figure size 400x400 with 0 Axes>
<Figure size 2000x2000 with 0 Axes>
df_nerve = pd.read_csv(f'{root}/results/{image_id}/data_nerve.csv', index_col=0)
df_nerve
| id | is_wsi | x_px | y_px | scale | total_fas | n_fibers | total_area | total_density | std_density | |
|---|---|---|---|---|---|---|---|---|---|---|
| 0 | sample_wsi_scale_0273 | True | 7813 | 7863 | 0.273 | 10 | 4155 | 0.583923 | 7115.665869 | 1005.693849 |
df_fas = pd.read_csv(f'{root}/results/{image_id}/data_fas.csv', index_col=0)
df_fas
| n_fas | n_fib | area | density | |
|---|---|---|---|---|
| 0 | 0 | 1119 | 0.182891 | 6118.394700 |
| 1 | 1 | 383 | 0.050129 | 7640.306931 |
| 2 | 2 | 946 | 0.140433 | 6736.305686 |
| 3 | 3 | 534 | 0.065160 | 8195.160650 |
| 4 | 4 | 366 | 0.043105 | 8490.945248 |
| 5 | 5 | 125 | 0.015002 | 8332.047390 |
| 6 | 6 | 164 | 0.023515 | 6974.354989 |
| 7 | 7 | 337 | 0.036705 | 9181.195894 |
| 8 | 8 | 106 | 0.016642 | 6369.273575 |
| 9 | 9 | 75 | 0.010340 | 7253.538672 |
df_fib = pd.read_csv(f'{root}/results/{image_id}/data_fib.csv', index_col=0)
df_fib = df_fib.rename(columns={'diameter_out': 'diameter'})
df_fib
| area_out | area_in | perimeter | circularity | convexity | solidity | eccentricity | major_axis_length | minor_axis_length | angle | diameter | diameter_in | thickness | g_ratio | n_fas | xmin | ymin | xmax | ymax | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 7.677978 | 3.292691 | 231.036577 | 0.606331 | 0.933163 | 0.906069 | 2.445910 | 5.186807 | 2.120604 | 34.596996 | 3.126644 | 2.047531 | 0.539556 | 0.654865 | 0.0 | 880.977778 | 205.297297 | 901.511111 | 227.675676 |
| 1 | 118.794754 | 39.889411 | 830.482315 | 0.726041 | 0.897158 | 0.973624 | 1.431737 | 14.793043 | 10.332236 | 97.859131 | 12.298544 | 7.126624 | 2.585960 | 0.579469 | 0.0 | 892.711111 | 172.216216 | 953.333333 | 218.918919 |
| 2 | 18.742553 | 9.934716 | 370.676186 | 0.574995 | 0.876249 | 0.875871 | 1.958578 | 7.144187 | 3.647640 | 174.681580 | 4.885055 | 3.556582 | 0.664237 | 0.728054 | 0.0 | 983.644444 | 191.675676 | 1006.133333 | 223.783784 |
| 3 | 105.212589 | 33.983733 | 748.825463 | 0.790917 | 0.916607 | 0.979055 | 1.222804 | 12.935622 | 10.578653 | 123.923126 | 11.574145 | 6.577951 | 2.498097 | 0.568331 | 0.0 | 929.866667 | 141.081081 | 980.711111 | 188.756757 |
| 4 | 27.972224 | 15.121934 | 386.475176 | 0.789420 | 0.921415 | 0.967320 | 1.252309 | 6.724402 | 5.369602 | 130.612839 | 5.967859 | 4.387920 | 0.789969 | 0.735259 | 0.0 | 1129.333333 | 154.702703 | 1159.644444 | 184.864865 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 4150 | 32.011696 | 8.365135 | 413.161468 | 0.790485 | 0.932958 | 0.977115 | 1.480821 | 7.836022 | 5.291674 | 131.498795 | 6.384243 | 3.263560 | 1.560342 | 0.511190 | 9.0 | 217.846154 | 320.307692 | 250.153846 | 351.692308 |
| 4151 | 46.576153 | 14.500362 | 502.173661 | 0.778539 | 0.932473 | 0.970675 | 1.495623 | 9.482864 | 6.340413 | 105.341042 | 7.700818 | 4.296794 | 1.702012 | 0.557966 | 9.0 | 422.769231 | 267.692308 | 463.384615 | 298.153846 |
| 4152 | 47.135121 | 17.149123 | 508.960457 | 0.767010 | 0.927581 | 0.956706 | 1.400718 | 9.372693 | 6.691348 | 97.500267 | 7.746890 | 4.672787 | 1.537051 | 0.603182 | 9.0 | 424.615385 | 237.230769 | 463.384615 | 268.615385 |
| 4153 | 4.652100 | 1.997377 | 169.095453 | 0.685819 | 0.920109 | 0.919022 | 1.585113 | 3.185095 | 2.009380 | 93.401825 | 2.433770 | 1.594722 | 0.419524 | 0.655248 | 9.0 | 430.153846 | 230.769231 | 448.615385 | 244.615385 |
| 4154 | 9.256502 | 4.119963 | 262.509666 | 0.566214 | 0.899516 | 0.936369 | 2.530454 | 5.572750 | 2.202273 | 88.055771 | 3.433037 | 2.290349 | 0.571344 | 0.667150 | 9.0 | 257.538462 | 390.461538 | 283.384615 | 408.000000 |
4155 rows × 19 columns
img = Image.open(f'{save_dir}/masked_fascicles.jpg')
plt.imshow(img, aspect='auto')
plt.axis('off')
(-0.5, 7862.5, 7812.5, -0.5)
col = 3
row = len(df_fas)//col + 1
plt.figure(figsize=(col*3, row*3))
plt.subplots_adjust(wspace=0.1, hspace=0.1)
for num in range(len(df_fas)):
plt.subplot(row, col, num+1)
img = Image.open(f'{save_dir}/fascicles/fascicle_{num:03}_bbox.jpg')
plt.imshow(img)
# plt.title(f'fascicle {num}', fontsize=10)
plt.axis('off')
fas_to_remove = [] # fill number(s) of improperly detected fas in the list
df_fib_rm = df_fib[~df_fib['n_fas'].isin(fas_to_remove)].dropna(axis=0)
df_fib_rm
| area_out | area_in | perimeter | circularity | convexity | solidity | eccentricity | major_axis_length | minor_axis_length | angle | diameter | diameter_in | thickness | g_ratio | n_fas | xmin | ymin | xmax | ymax | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 7.677978 | 3.292691 | 231.036577 | 0.606331 | 0.933163 | 0.906069 | 2.445910 | 5.186807 | 2.120604 | 34.596996 | 3.126644 | 2.047531 | 0.539556 | 0.654865 | 0.0 | 880.977778 | 205.297297 | 901.511111 | 227.675676 |
| 1 | 118.794754 | 39.889411 | 830.482315 | 0.726041 | 0.897158 | 0.973624 | 1.431737 | 14.793043 | 10.332236 | 97.859131 | 12.298544 | 7.126624 | 2.585960 | 0.579469 | 0.0 | 892.711111 | 172.216216 | 953.333333 | 218.918919 |
| 2 | 18.742553 | 9.934716 | 370.676186 | 0.574995 | 0.876249 | 0.875871 | 1.958578 | 7.144187 | 3.647640 | 174.681580 | 4.885055 | 3.556582 | 0.664237 | 0.728054 | 0.0 | 983.644444 | 191.675676 | 1006.133333 | 223.783784 |
| 3 | 105.212589 | 33.983733 | 748.825463 | 0.790917 | 0.916607 | 0.979055 | 1.222804 | 12.935622 | 10.578653 | 123.923126 | 11.574145 | 6.577951 | 2.498097 | 0.568331 | 0.0 | 929.866667 | 141.081081 | 980.711111 | 188.756757 |
| 4 | 27.972224 | 15.121934 | 386.475176 | 0.789420 | 0.921415 | 0.967320 | 1.252309 | 6.724402 | 5.369602 | 130.612839 | 5.967859 | 4.387920 | 0.789969 | 0.735259 | 0.0 | 1129.333333 | 154.702703 | 1159.644444 | 184.864865 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 4150 | 32.011696 | 8.365135 | 413.161468 | 0.790485 | 0.932958 | 0.977115 | 1.480821 | 7.836022 | 5.291674 | 131.498795 | 6.384243 | 3.263560 | 1.560342 | 0.511190 | 9.0 | 217.846154 | 320.307692 | 250.153846 | 351.692308 |
| 4151 | 46.576153 | 14.500362 | 502.173661 | 0.778539 | 0.932473 | 0.970675 | 1.495623 | 9.482864 | 6.340413 | 105.341042 | 7.700818 | 4.296794 | 1.702012 | 0.557966 | 9.0 | 422.769231 | 267.692308 | 463.384615 | 298.153846 |
| 4152 | 47.135121 | 17.149123 | 508.960457 | 0.767010 | 0.927581 | 0.956706 | 1.400718 | 9.372693 | 6.691348 | 97.500267 | 7.746890 | 4.672787 | 1.537051 | 0.603182 | 9.0 | 424.615385 | 237.230769 | 463.384615 | 268.615385 |
| 4153 | 4.652100 | 1.997377 | 169.095453 | 0.685819 | 0.920109 | 0.919022 | 1.585113 | 3.185095 | 2.009380 | 93.401825 | 2.433770 | 1.594722 | 0.419524 | 0.655248 | 9.0 | 430.153846 | 230.769231 | 448.615385 | 244.615385 |
| 4154 | 9.256502 | 4.119963 | 262.509666 | 0.566214 | 0.899516 | 0.936369 | 2.530454 | 5.572750 | 2.202273 | 88.055771 | 3.433037 | 2.290349 | 0.571344 | 0.667150 | 9.0 | 257.538462 | 390.461538 | 283.384615 | 408.000000 |
3893 rows × 19 columns
len(df_fib_rm[df_fib_rm['diameter'] >= 25]), df_fib_rm.isna().sum().sum()
(0, 0)
df_fas_rm = df_fas[~df_fas['n_fas'].isin(fas_to_remove)]
df_fas_rm
| n_fas | n_fib | area | density | |
|---|---|---|---|---|
| 0 | 0 | 1119 | 0.182891 | 6118.394700 |
| 1 | 1 | 383 | 0.050129 | 7640.306931 |
| 2 | 2 | 946 | 0.140433 | 6736.305686 |
| 3 | 3 | 534 | 0.065160 | 8195.160650 |
| 4 | 4 | 366 | 0.043105 | 8490.945248 |
| 5 | 5 | 125 | 0.015002 | 8332.047390 |
| 6 | 6 | 164 | 0.023515 | 6974.354989 |
| 7 | 7 | 337 | 0.036705 | 9181.195894 |
| 8 | 8 | 106 | 0.016642 | 6369.273575 |
| 9 | 9 | 75 | 0.010340 | 7253.538672 |
df_nerve_rm = df_nerve.copy()
df_nerve_rm['total_fas'] = len(df_fas_rm)
df_nerve_rm['n_fibers'] = df_fas_rm['n_fib'].sum()
df_nerve_rm['total_area'] = df_fas_rm['area'].sum()
df_nerve_rm['total_density'] = df_nerve_rm['n_fibers'] / df_nerve_rm['total_area']
df_nerve_rm['std_density'] = df_fas_rm['density'].std()
df_nerve_rm
| id | is_wsi | x_px | y_px | scale | total_fas | n_fibers | total_area | total_density | std_density | |
|---|---|---|---|---|---|---|---|---|---|---|
| 0 | sample_wsi_scale_0273 | True | 7813 | 7863 | 0.273 | 10 | 4155 | 0.583923 | 7115.665869 | 1005.693849 |
def get_density(df, area, name, xlim=20, dev=1):
idx = []
density = []
bins = int(xlim // dev) + 1
for i in range(bins):
if i == bins - 1:
dens = len(df[df[name] >= (i) * dev]) / area
else:
dens = len(df[ (df[name] >= (i) * dev) & (df[name] < (i + 1) * dev) ]) / area
idx.append(i * dev)
density.append(dens)
density_df = pd.DataFrame({name: density}, index = idx)
return density_df
def draw_histogram(df, total_area, name='diameter', dev=1.0, width=1, xlim=20, tick=1):
density = get_density(df, total_area, name=name, xlim=xlim, dev=dev)
fig, ax = plt.subplots(figsize=(6, 3))
# ax.set_title(name)
ax.set_xlim(0, xlim)
ax.set_xlabel(f'{name} [um]', fontsize=10)
ax.set_ylabel('Density [/mm2]', fontsize=10)
ax.set_xticks(np.arange(0, xlim, tick))
ax.set_xticks(np.arange(0, xlim, dev), minor=True)
ax.bar(density.index, density[name], width=width, color='teal')
# ax.bar(density_ref.index+width/2, density_ref[name], width=width, alpha=0.5, color='powderblue')
# ax.legend(['sample', 'normal 57 y.o. male'])
plt.savefig(f'{save_dir}/histogram_{name}.jpg')
area = df_nerve_rm['total_area'][0]
draw_histogram(df_fib_rm, area, name='diameter', dev=1, width=1, xlim=20, tick=1)
draw_histogram(df_fib_rm, area, name='thickness', dev=0.20, width=0.20, xlim=5, tick=1)
palette = sns.color_palette('gray', n_colors=2)
sns.set(context='talk', style='ticks', palette=palette, font_scale=.7,rc={
'axes.spines.right': False,
'axes.spines.top': False,
'figure.figsize': [4.0, 4.0],
'legend.fancybox': False,
'legend.edgecolor': 'black',
'legend.fancybox': False,
'legend.frameon': True,
'legend.handletextpad': 0.5,
'legend.loc': 'upper right',
'legend.loc': 'best',
'legend.title_fontsize': None,
})
sns.palplot(palette)
sns.histplot(df_fib_rm['diameter'], bins=40)
plt.xlim(0, 20)
plt.xlabel('Diameter [um]', fontsize=10)
plt.ylabel('Count / slide', fontsize=10)
plt.savefig(f'{save_dir}/histogram_diameter.jpg')
sns.histplot(df_fib_rm['thickness'], bins=50)
plt.xlim(0, 5)
plt.xlabel('Thickness [um]', fontsize=10)
plt.ylabel('Count / slide', fontsize=10)
plt.savefig(f'{save_dir}/histogram_thickness.jpg')
sns.scatterplot(data=df_fib_rm, x='diameter', y='thickness', s=2)
plt.xlabel('Diameter [um]', fontsize=10)
plt.ylabel('Thickness [um]', fontsize=10)
plt.xlim(0, 20)
plt.ylim(0, 5)
plt.savefig(f'{save_dir}/scatter_diameter_thickness.jpg')
sns.scatterplot(data=df_fib_rm, x='diameter', y='g_ratio', s=2)
plt.xlim(0, 20)
plt.ylim(0, 1)
plt.savefig(f'{save_dir}/scatter_diameter_g_ratio.jpg')
gmm = GaussianMixture(n_components=2, covariance_type='tied')
loaded_gmm = joblib.load(f'{root}/weights/gmm_classifier.pkl')
df_fib_rm['gmm'] = loaded_gmm.predict(df_fib_rm[['diameter', 'thickness']])
palette = sns.color_palette('viridis', n_colors=3)
sns.set(context='talk', style='ticks', palette=palette, font_scale=.7,rc={
'axes.spines.right': False,
'axes.spines.top': False,
'figure.figsize': [4.0, 4.0],
'legend.fancybox': False,
'legend.edgecolor': 'black',
'legend.fancybox': False,
'legend.frameon': True,
'legend.handletextpad': 0.5,
'legend.loc': 'upper right',
'legend.loc': 'best',
'legend.title_fontsize': None,
})
sns.palplot(palette)
sns.scatterplot(data=df_fib_rm, x='diameter', y='g_ratio', hue='gmm', s=2)
plt.legend(bbox_to_anchor=(1.05, 1), loc='upper left', borderaxespad=0.)
plt.xlabel('Diameter [um]', fontsize=10)
plt.ylabel('g-ratio', fontsize=10)
plt.xlim(0, 20)
plt.ylim(0, 1)
plt.savefig(f'{save_dir}/scatter_diameter_g_ratio_gmm.jpg')
sns.scatterplot(data=df_fib_rm, x='diameter', y='thickness', hue='gmm', s=2)
plt.legend(bbox_to_anchor=(1.05, 1), loc='upper left', borderaxespad=0.)
plt.xlabel('Diameter [um]', fontsize=10)
plt.ylabel('Thickness [um]', fontsize=10)
plt.xlim(0, 20)
plt.ylim(0, 5)
plt.savefig(f'{save_dir}/scatter_diameter_thickness_gmm.jpg')
sns.histplot(data=df_fib_rm, x='diameter', bins=40, hue='gmm')
labels = ['small fiber', 'large fiber']
plt.legend(labels=labels, bbox_to_anchor=(1.05, 1), loc='upper left', borderaxespad=0.)
plt.xlim(0, 20)
plt.xlabel('Diameter [um]', fontsize=10)
plt.ylabel('Count / slide', fontsize=10)
plt.savefig(f'{save_dir}/histogram_diameter_gmm.jpg')
sns.histplot(data=df_fib_rm, x='thickness', bins=50, hue='gmm')
labels = ['small fiber', 'large fiber']
plt.legend(labels=labels, bbox_to_anchor=(1.05, 1), loc='upper left', borderaxespad=0.)
plt.xlim(0, 5)
plt.xlabel('Thickness [um]', fontsize=10)
plt.ylabel('Count / slide', fontsize=10)
plt.savefig(f'{save_dir}/histogram_thickness_gmm.jpg')
def get_linear_regression(df, xlabel='diameter', ylabel='g_ratio'):
df = df.copy().dropna(subset=(xlabel, ylabel))
x = df[xlabel].values.reshape(-1, 1)
y = df[ylabel]
lr = LinearRegression()
lr.fit(x, y)
# r2 = lr.score(x, y)
coef = lr.coef_
intercept = lr.intercept_
sns.jointplot(data=df_fib_rm, x='diameter', xlim=(0, 15), ylim=(0, 1), y='g_ratio', hue='gmm', s=2)
plt.legend([],[], frameon=False)
x_ = np.arange(min(df[xlabel]), max(df[xlabel]), 0.1)
plt.plot(x_, coef*x_+intercept, color='gray', linewidth=2)
plt.xlim(0, 20)
plt.ylim(0, 1)
plt.xlabel('Diameter [um]', fontsize=10)
plt.ylabel('g-ratio', fontsize=10)
plt.savefig(f'{save_dir}/{xlabel}_{ylabel}_lr.jpg')
print(
f'coef:{coef[0]:.4f}, intercept:{intercept:.3f}'
)
return coef[0], intercept
df_nerve['coef'], df_nerve['intercept'] = get_linear_regression(df_fib_rm)
coef:-0.0156, intercept:0.723
df_fib
| area_out | area_in | perimeter | circularity | convexity | solidity | eccentricity | major_axis_length | minor_axis_length | angle | diameter | diameter_in | thickness | g_ratio | n_fas | xmin | ymin | xmax | ymax | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 7.677978 | 3.292691 | 231.036577 | 0.606331 | 0.933163 | 0.906069 | 2.445910 | 5.186807 | 2.120604 | 34.596996 | 3.126644 | 2.047531 | 0.539556 | 0.654865 | 0.0 | 880.977778 | 205.297297 | 901.511111 | 227.675676 |
| 1 | 118.794754 | 39.889411 | 830.482315 | 0.726041 | 0.897158 | 0.973624 | 1.431737 | 14.793043 | 10.332236 | 97.859131 | 12.298544 | 7.126624 | 2.585960 | 0.579469 | 0.0 | 892.711111 | 172.216216 | 953.333333 | 218.918919 |
| 2 | 18.742553 | 9.934716 | 370.676186 | 0.574995 | 0.876249 | 0.875871 | 1.958578 | 7.144187 | 3.647640 | 174.681580 | 4.885055 | 3.556582 | 0.664237 | 0.728054 | 0.0 | 983.644444 | 191.675676 | 1006.133333 | 223.783784 |
| 3 | 105.212589 | 33.983733 | 748.825463 | 0.790917 | 0.916607 | 0.979055 | 1.222804 | 12.935622 | 10.578653 | 123.923126 | 11.574145 | 6.577951 | 2.498097 | 0.568331 | 0.0 | 929.866667 | 141.081081 | 980.711111 | 188.756757 |
| 4 | 27.972224 | 15.121934 | 386.475176 | 0.789420 | 0.921415 | 0.967320 | 1.252309 | 6.724402 | 5.369602 | 130.612839 | 5.967859 | 4.387920 | 0.789969 | 0.735259 | 0.0 | 1129.333333 | 154.702703 | 1159.644444 | 184.864865 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 4150 | 32.011696 | 8.365135 | 413.161468 | 0.790485 | 0.932958 | 0.977115 | 1.480821 | 7.836022 | 5.291674 | 131.498795 | 6.384243 | 3.263560 | 1.560342 | 0.511190 | 9.0 | 217.846154 | 320.307692 | 250.153846 | 351.692308 |
| 4151 | 46.576153 | 14.500362 | 502.173661 | 0.778539 | 0.932473 | 0.970675 | 1.495623 | 9.482864 | 6.340413 | 105.341042 | 7.700818 | 4.296794 | 1.702012 | 0.557966 | 9.0 | 422.769231 | 267.692308 | 463.384615 | 298.153846 |
| 4152 | 47.135121 | 17.149123 | 508.960457 | 0.767010 | 0.927581 | 0.956706 | 1.400718 | 9.372693 | 6.691348 | 97.500267 | 7.746890 | 4.672787 | 1.537051 | 0.603182 | 9.0 | 424.615385 | 237.230769 | 463.384615 | 268.615385 |
| 4153 | 4.652100 | 1.997377 | 169.095453 | 0.685819 | 0.920109 | 0.919022 | 1.585113 | 3.185095 | 2.009380 | 93.401825 | 2.433770 | 1.594722 | 0.419524 | 0.655248 | 9.0 | 430.153846 | 230.769231 | 448.615385 | 244.615385 |
| 4154 | 9.256502 | 4.119963 | 262.509666 | 0.566214 | 0.899516 | 0.936369 | 2.530454 | 5.572750 | 2.202273 | 88.055771 | 3.433037 | 2.290349 | 0.571344 | 0.667150 | 9.0 | 257.538462 | 390.461538 | 283.384615 | 408.000000 |
4155 rows × 19 columns
# down scale to 1/10 to compute fast
df_fib['x'] = (df_fib['xmin'] + df_fib['xmax']) / 20
df_fib['y'] = (df_fib['ymin'] + df_fib['ymax']) / 20
df_fas_rm['NND'] = np.nan
df_fas_rm['KDE_entropy'] = np.nan
for idx in range(len(df_fas_rm)):
data = df_fib[(df_fib['n_fas']==idx)].dropna(subset=['x', 'y'])
data = data[['x', 'y']]
x = data['x']
y = data['y']
if len(data) < 10:
continue
# Calculate nearest neighbor distances
distances = distance.cdist(data, data)
nnd = np.min(np.ma.masked_array(distances, mask=np.eye(len(data))), axis=1)
# Find the minimum distance for each point
masked_distances = np.ma.masked_array(distances, mask=np.eye(len(data)))
nearest_distances = masked_distances.min(axis=1)
# Calculate the mean nearest neighbor distance
mean_nnd = nearest_distances.mean() * 10 * scale # um
df_fas_rm['NND'][idx] = mean_nnd
print(f"Mean Nearest Neighbor Distance: {mean_nnd:.2f} um")
# Perform kernel density estimation
bandwidth = 5 # note that 1/10 scale to original pixel
kde = KernelDensity(bandwidth=bandwidth, kernel='gaussian')
kde.fit(data)
# Evaluate the KDE on the grid
grid_size = 1
grid_x, grid_y = np.mgrid[0:x.max():grid_size, 0:y.max():grid_size]
grid_points = np.column_stack((grid_x.ravel(), grid_y.ravel()))
log_density = kde.score_samples(grid_points)
density = np.exp(log_density)
# Calculate the entropy of the KDE distribution
entropy_value = entropy(density)
df_fas_rm['KDE_entropy'][idx] = entropy_value
print(f"Entropy: {entropy_value:.2f}")
# plot data
fig, ax = plt.subplots(1, 4, figsize=(12, 3))
img = Image.open(f'{save_dir}/fascicles/fascicle_{idx:03}_bbox.jpg')
ax[0].imshow(img)
ax[1].set_aspect('equal')
ax[1].scatter(x, y, s=5)
ax[1].set_xlim(0, x.max()+10)
ax[1].set_ylim(0, y.max()+10)
ax[1].invert_yaxis()
ax[1].axis('off')
ax[1].set_aspect('equal')
im = ax[2].imshow(density.reshape(grid_x.shape).T, origin='lower', cmap='viridis')
ax[2].invert_yaxis()
ax[2].axis('off')
ax[2].set_aspect('equal')
# colorbar
cbar = fig.colorbar(im, ax=ax[3], pad=0.05, location='left')
cbar.set_ticks([])
ax[3].axis('off')
plt.savefig(f'{save_dir}/fascicles/fascicle_{idx:03}_kde.jpg')
# plt.close()
df_fas_rm['expected_nnd'] = 1 / (2 * np.sqrt(df_fas_rm['n_fib'] / (df_fas_rm['area'] * 1000000))) # um
# Calculate the ratio of observed mean nearest neighbor distance to expected value
df_fas_rm['NNI'] = df_fas_rm['NND'] / df_fas_rm['expected_nnd']
Mean Nearest Neighbor Distance: 8.13 um Entropy: 10.10 Mean Nearest Neighbor Distance: 7.50 um Entropy: 8.89 Mean Nearest Neighbor Distance: 8.12 um Entropy: 9.87 Mean Nearest Neighbor Distance: 7.04 um Entropy: 9.13 Mean Nearest Neighbor Distance: 7.16 um Entropy: 8.75 Mean Nearest Neighbor Distance: 6.81 um Entropy: 7.64 Mean Nearest Neighbor Distance: 7.88 um Entropy: 8.12 Mean Nearest Neighbor Distance: 6.77 um Entropy: 8.59 Mean Nearest Neighbor Distance: 8.01 um Entropy: 7.86 Mean Nearest Neighbor Distance: 7.07 um Entropy: 7.38
idx = 0
# for idx in range(len(df_fas_rm)):
data = df_fib[(df_fib['n_fas']==idx)].dropna(subset=['x', 'y'])
data = data[['x', 'y']]
x = data['x']
y = data['y']
# if len(data) < 10:
# continue
# Calculate nearest neighbor distances
distances = distance.cdist(data, data)
nnd = np.min(np.ma.masked_array(distances, mask=np.eye(len(data))), axis=1)
# Find the minimum distance for each point
masked_distances = np.ma.masked_array(distances, mask=np.eye(len(data)))
nearest_distances = masked_distances.min(axis=1)
# Calculate the mean nearest neighbor distance
mean_nnd = nearest_distances.mean() * 10 * scale # um
df_fas_rm['NND'][idx] = mean_nnd
print(f"Mean Nearest Neighbor Distance: {mean_nnd:.2f} um")
# Perform kernel density estimation
bandwidth = 5 # note that 1/10 scale to original pixel
kde = KernelDensity(bandwidth=bandwidth, kernel='gaussian')
kde.fit(data)
# Evaluate the KDE on the grid
grid_size = 1
grid_x, grid_y = np.mgrid[0:x.max():grid_size, 0:y.max():grid_size]
grid_points = np.column_stack((grid_x.ravel(), grid_y.ravel()))
log_density = kde.score_samples(grid_points)
density = np.exp(log_density)
# Calculate the entropy of the KDE distribution
entropy_value = entropy(density)
df_fas_rm['KDE_entropy'][idx] = entropy_value
print(f"Entropy: {entropy_value:.2f}")
# plot data
fig, ax = plt.subplots(1, 4, figsize=(12, 3))
img = Image.open(f'{save_dir}/fascicles/fascicle_{idx:03}_bbox.jpg')
ax[0].imshow(img)
ax[1].set_aspect('equal')
ax[1].scatter(x, y, s=5)
ax[1].set_xlim(0, x.max()+10)
ax[1].set_ylim(0, y.max()+10)
ax[1].invert_yaxis()
ax[1].axis('off')
ax[1].set_aspect('equal')
im = ax[2].imshow(density.reshape(grid_x.shape).T, origin='lower', cmap='viridis')
ax[2].invert_yaxis()
ax[2].axis('off')
ax[2].set_aspect('equal')
# colorbar
cbar = fig.colorbar(im, ax=ax[3], pad=0.05, location='left')
cbar.set_ticks([])
ax[3].axis('off')
plt.savefig(f'{save_dir}/fascicles/fascicle_{idx:03}_kde.jpg')
# plt.close()
Mean Nearest Neighbor Distance: 8.13 um Entropy: 10.10
density.min(), density.max(), density.mean()
(4.901244531291384e-31, 7.50410256471614e-05, 2.9131704750970153e-05)
density.sum()
0.987790333117873
entropy(density[:4])
0.5612423377111373
density.shape
(1800,)
df_fas_rm
| n_fas | n_fib | area | density | NND | KDE_entropy | expected_nnd | NNI | |
|---|---|---|---|---|---|---|---|---|
| 0 | 0 | 1119 | 0.182891 | 6118.394700 | 8.127743 | 10.104460 | 6.392213 | 1.271507 |
| 1 | 1 | 383 | 0.050129 | 7640.306931 | 7.495048 | 8.894850 | 5.720245 | 1.310267 |
| 2 | 2 | 946 | 0.140433 | 6736.305686 | 8.118978 | 9.869836 | 6.091989 | 1.332730 |
| 3 | 3 | 534 | 0.065160 | 8195.160650 | 7.043000 | 9.134708 | 5.523206 | 1.275165 |
| 4 | 4 | 366 | 0.043105 | 8490.945248 | 7.158880 | 8.746445 | 5.426152 | 1.319329 |
| 5 | 5 | 125 | 0.015002 | 8332.047390 | 6.814293 | 7.636691 | 5.477648 | 1.244018 |
| 6 | 6 | 164 | 0.023515 | 6974.354989 | 7.878395 | 8.122971 | 5.987120 | 1.315891 |
| 7 | 7 | 337 | 0.036705 | 9181.195894 | 6.773026 | 8.586368 | 5.218196 | 1.297963 |
| 8 | 8 | 106 | 0.016642 | 6369.273575 | 8.008600 | 7.862696 | 6.265057 | 1.278296 |
| 9 | 9 | 75 | 0.010340 | 7253.538672 | 7.070057 | 7.382267 | 5.870770 | 1.204281 |
df_nerve_rm
| id | is_wsi | x_px | y_px | scale | total_fas | n_fibers | total_area | total_density | std_density | |
|---|---|---|---|---|---|---|---|---|---|---|
| 0 | sample_wsi_scale_0273 | True | 7813 | 7863 | 0.273 | 10 | 4155 | 0.583923 | 7115.665869 | 1005.693849 |
df_nerve_rm.to_csv(f'{save_dir}/df_nerve_rm.csv')
df_fas_rm.to_csv(f'{save_dir}/df_fas_rm.csv')
df_fib_rm.to_csv(f'{save_dir}/df_fib_rm.csv')
# from google.colab import runtime
# runtime.unassign()